#include <cxcore.h>
#include <cv.h>
#include <highgui.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#define PI 3.14159265
int main(){
    
    IplImage * resim;
    CvScalar s,s1;
    int i,j,Y,Cb,Cr;
    resim = cvLoadImage("lc.jpg",1);
    IplImage* ycb = cvCreateImage(cvSize(resim->width,resim->height), IPL_DEPTH_8U, 3);
    IplImage* skin = cvCreateImage(cvSize(resim->width,resim->height), IPL_DEPTH_8U, 3);
    IplImage* mask = cvCreateImage(cvSize(resim->width,resim->height), IPL_DEPTH_8U, 1);
    cvCvtColor(resim,ycb,CV_BGR2YCrCb);
    
    CvConnectedComp comp;
    CvScalar floodColor;
    CvPoint floodSeed;
    int lo=1;
    int up=1;
    floodColor = cvRealScalar(128);
    floodSeed=cvPoint(10,10);
    
    CvMemStorage * storage = cvCreateMemStorage(0);
    CvSeq * contour = 0;
    int mode = CV_RETR_EXTERNAL;
    mode = CV_RETR_CCOMP;
    
    

    CvMat *cb = cvCreateMat( resim->height, resim->width, CV_32FC1);
    CvMat *cr = cvCreateMat( resim->height, resim->width, CV_32FC1);
    
    const unsigned char Kl = 125;
    const unsigned char Kh = 188;
    const float Wcb = 46.97;
    const unsigned char WLcb = 23;
    const unsigned char WHcb = 14;
    const float Wcr = 38.76;
    const unsigned char WLcr = 20;
    const unsigned char WHcr = 10;
    const unsigned char Ymin = 16;
    const unsigned char Ymax = 235;    
    float CbP, CrP,WcbY,WcrY,CbBarY,CrBarY;
    float a2=644.6521;
    float b2=196.8409,x,y;
   for(i=0;i<(ycb->height)-1;i++){
      for(j=0;j<(ycb->width)-1;j++){
         s=cvGet2D(ycb,i,j);
         Y=(int)s.val[0];
         Cr=(int)s.val[1];
         Cb=(int)s.val[2];
         if (Y >= Kl && Y <= Kh){
          CbP = Cb;
          CrP = Cr;
         }else{
             if (Y < Kl){
              WcbY = WLcb + ((Y - Ymin) * (Wcb - WLcb)) / (Kl - Ymin);
              WcrY = WLcr + ((Y - Ymin) * (Wcr - WLcr)) / (Kl - Ymin);
    
              CbBarY = 108 + ((Kl - Y) * (118 - 108)) / (Kl - Ymin);
              CrBarY = 154 - ((Kl - Y) * (154 - 144)) / (Kl - Ymin);
             }else{ 
              WcbY = WHcb + ((Ymax - Y) * (Wcb - WHcb)) / (Ymax - Kh);
              WcrY = WHcr + ((Ymax - Y) * (Wcr - WHcr)) / (Ymax - Kh);

              CbBarY = 108 + ((Y - Kh) * (118 - 108)) / (Ymax - Kh);
              CrBarY = 154 - ((Y - Kh) * (154 - 132)) / (Ymax - Kh);
             }
              CbP = (Cb - CbBarY) * (Wcb / WcbY) + CbBarY;
              CrP = (Cr - CrBarY) * (Wcr / WcrY) + CrBarY;
          }
         
    //      cvmSet(cb,i,j,CbP);
    //      cvmSet(cr,i,j,CrP);
          s.val[2]=(CbP < 0) ? 0 : ((CbP > 255) ? 255: static_cast <BYTE> (CbP));
          s.val[1]=(CrP < 0) ? 0 : ((CrP > 255) ? 255: static_cast <BYTE> (CrP));
          CbP=(float)s.val[2];
          CrP=(float)s.val[1];
          
          x=(cos(2.53)*(CbP-109.38)+ sin(2.53)*(CrP-152.02));
          y=((-1)*sin(2.53)*(CbP-109.38)+ cos(2.53)*(CrP-152.02));
          //printf("%d %d\n",(int)CbP,(int)CrP);
         // fprintf(fp,"%f %f\n",x,y);
          if(((pow((x-1.6),2)/a2)+(pow((y-2.41),2)/b2))<=1.0){
          s1.val[0]=255;                                               
          }else{
          s1.val[0]=0;
          } 

          cvSet2D(skin,i,j,s);
          cvSet2D(mask,i,j,s1);
          
          
       }
    }
    cvNamedWindow( "result1", 1 );
    IplConvKernel* str = cvCreateStructuringElementEx(3, 3, 1, 1, CV_SHAPE_ELLIPSE, NULL) ;
  //  cvDilate(mask,mask,NULL,2);
  //  cvErode(mask,mask,NULL,1);
            cvShowImage( "result1",mask );


    cvMorphologyEx(mask,mask,NULL,str,CV_MOP_CLOSE,3);
   // cvMorphologyEx(mask,mask,NULL,NULL,CV_MOP_CLOSE,5);

    //cvFindContours(mask, storage, &contour, sizeof(CvContour), mode, CV_CHAIN_APPROX_SIMPLE, cvPoint(0,0));
//draw the contour
   // cvDrawContours(mask, contour, cvRealScalar(100), cvRealScalar(179), 2, 2, 8);
   // cvFloodFill( mask, floodSeed, floodColor, cvRealScalar(lo), cvRealScalar(up), &comp, 8, NULL);
    cvNamedWindow( "result", 1 );
    cvShowImage( "result",mask );
    cvWaitKey( -1 );
    cvDestroyWindow("result");
    cvDestroyWindow("result1");
    cvReleaseMat(&cb);
    cvReleaseMat(&cr);

    cvSaveImage("mask.jpg",mask);
    cvSaveImage("skin.jpg",skin);
        cvReleaseImage( &mask );
            cvReleaseImage( &skin );
                cvReleaseImage( &ycb);
    cvReleaseStructuringElement(&str);
    system("pause");
    return 0;
    }
